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Abstract. We study the problem of a quantum quench in which the initial state is 
the ground state of an inhomogeneous hamiltonian, in two different models, conformal 
field theory and ordinary free field theory, which are known to exhibit thermalisation of 
finite regions in the homogeneous case. We derive general expressions for the evolution 
of the energy flow and correlation functions, as well as the entanglement entropy in 
the conformal case. Comparison of the results of the two approaches in the regime of 
their common validity shows agreement up to a point further discussed. Unlike the 
thermal analogue, the evolution in our problem is non-diffusive and can be physically 
interpreted using an intuitive picture of quasiparticlcs emitted from the initial time 
hypersurface and propagating semiclassically. 
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1. Introduction 

An isolated quantum system in which some of the parameters that determine the 
dynamics change rapidly at a specific instant is said to undergo a quantum quench. 
Recently such rapid changes have become experimentally feasible in cold atom systems 
[1, 2, 3, 4]. These experiments as well as the development of numerical techniques 
for the study of quantum dynamics (t-DMRG) have motivated significant theoretical 
interest as the investigation of a wide range of problems, either in lattice models 
[5, 6, 7, 8, 9, 10, 11, 12, 13, 14] or in continuous theories [7, 8, 15] shows. The 
importance of quantum quenches relies on the potential discovery of novel physical 
phenomena and the fact that they will help us to better understand out of equilibrium 
quantum behaviour. On the other hand they can be described by simple theoretical 
models. More specifically, to find the evolution that follows a quantum quench in the 
Heisenberg picture one should calculate the expectation values of operators whose time 
dependence is determined by the hamiltonian after the quench, with respect to the pure 
state in which the system was before the quench, which is typically the ground state 
of the initial hamiltonian. Then one finds that while in systems with finite number of 
degrees of freedom local observables exhibit periodic or quasiperiodic behaviour, in the 
thermodynamic limit they tend to stationary values instead. It should be emphasised 
that this happens to connected correlation functions of local observables, as a result of 
the interference of the infinite number of momentum modes of the whole system. 

In particular let us consider a system of coupled harmonic oscillators or equivalently 
a free field theory, described by a general dispersion relation with some energy gap or 
"mass" m and maximum group velocity of excitations v max . Assume that the system 
lies on the ground state of the initial hamiltonian H when at time t — the mass 
is quenched from m to a different value m ^ 0. Then after the quench there is an 
extensive excess in energy in comparison with the ground state of the final hamiltonian 
H, which is distributed to the excitation levels of H . In a spacetime representation 
these excitations appear as quasiparticles that emerge from the t — hypersurface 
and propagate forward in time. For free field theories the two-point correlation function 
contains all the information required to determine their state since this is a superposition 
of gaussian wavepackets. This correlation function of two points separated by distance r 
turns out [7, 8] to remain unaffected by the quench until time t = r/2v max when it starts 
changing and finally for large times it takes the form of a thermal correlation function 
with a momentum dependent effective temperature. Additionally if mo > m then for 
large values of m the effective temperature l//3 e // is of order m and asymptotically 
independent of the momentum [16] . To understand this behaviour we can imagine pairs 
of coherent quasiparticles emitted by neighbouring points of the t — hypersurface that 
induce correlations between spatially separated points as soon as the fastest ones reach 
them after time t = r/2v max - This is called the horizon effect [7, 8]. On the other hand 
for large times, the interference between a large number of incoherent quasiparticles 
coming from different and uncorrelated points of the t — hypersurface drives the 
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system to stationary and in particular thermal behaviour as far as local observables are 
concerned. This means that any finite subsystem tends to thermal equilibrium with 
its complement which acts as a bath {thermalisation) [17, 18, 7, 8, 12] (also [19, 20] 
for recent general discussions of the subject). The fact that the effective temperature 
depends on the momentum k in the free field theory case should be expected since 
the final hamiltonian can be diagonalised exactly in momentum space which means 
that the different momentum modes evolve independently and so thermalise to different 
temperatures. 

The above observations have been shown to be valid not only in the simple case of a 
free field theory but also in two other important general cases. First [7, 8] in conformal 
field theory (CFT), which describes one-dimensional quantum systems at criticality 
in the continuum limit (equivalently massless interacting field theories). Second [21] 
in a self-consistent Hartree-Fock approximation of a system of anharmonic oscillators 
(equivalently an interacting bosonic field theory), which is valid in the large- N limit of 
the linear sigma model. 

An important question is to what extent the stationary behaviour of the system 
depends on the initial state and if there are any physical quantities other than the 
conserved energy (and possibly any other quantity that commutes with the hamiltonian) 
about which information survives in the final state. One way to study the effect of 
different initial states on the evolution is to choose the hamiltonian before the quench 
to be spatially inhomogeneous. We call this kind of quench an inhomogeneous quench. 
More specifically we will concentrate on the previous system of harmonic oscillators 
where the initial mass mo now depends on the position x. In this case the initial 
hamiltonian is not diagonal in momentum space although it is still diagonalizable. A 
question of particular interest is whether the evolution after the quench resembles that 
of the thermal analogue, which would be the heat diffusion occurring in a system with 
inhomogeneous initial temperature distribution u (x). A special case of distribution 
where this comparison should be easier is the step distribution, that is when uq(x), 
respectively m (x), has different values for x > and x < 0. In the thermal case the 
diffusive nature of the heat equation leads to the heat current across x = decreasing 
with time as l/y/i and we wish to study what happens following an equivalent initial 
state in the case of a quantum quench. To this end we will be calculating the energy 
flow instead, since there is not such a notion as the heat current in our case. Recall 
that in the thermal problem the heat current j q is related to the energy current j e and 
the particle current j n according to j q = j e — /ij n where \i is the chemical potential, and 
that these currents also exhibit the same diffusive behaviour. 

At the same time we will derive general expressions for the correlation functions 
and the entanglement entropy which measure the correlations between points or parts 
of the system respectively. Similar calculations in the special case of a step distribution 
and m = have been done using a CFT method alternative to ours in [22]. Also local 
effects caused by defects have been investigated in [9] for a general lattice model and in 
[15] for a continuous field theory. General analytical results for the calculation of the 
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entanglement entropy in field theory can be found in [23, 24, 25, 26]. 

In this paper we consider a continuous bosonic field with a relativistic dispersion 
relation to = \Jk 2 + m? where v max = 1, although other dispersion relations drawn 
from lattice models and having the same energy gap and maximum group velocity 
are expected to lead to similar behaviour. For our purposes we can choose the initial 
distribution to depend only on one space coordinate and so it is sufficient to consider 
only one dimensional systems. Also we often make use of the so-called deep quench 
limit, that is the limit when m m. This must obviously reflect all the characteristic 
features of a quantum quench since it is one of the two most extreme possibilities for 
the relation between the two masses. In our inhomogeneous problem this means that 
m should be much smaller than any value of the initial mass distribution itlq{x). 

This paper is organised as follows: in section 2 we apply the methods of CFT to 
solve the general problem in the case where the theory after the quench is massless. In 
section 3 we analyse two special cases where the initial distribution has a bump or a 
step. In section 4 we use free field theory methods to solve the massive problem as well 
and in section 5 we compare the results obtained from these two approaches. Lastly 
in section 6 we discuss our findings as compared to the thermal analogue and give a 
physical interpretation. 

2. Massless case - CFT approach 

As shown in earlier work [8], the problem of a quantum quench can be mapped to a 
Euclidean field theory defined on a strip (or a rf+l-dimensional slab in general) where 
the transverse direction corresponds to imaginary time r. To see this let us consider the 
expectation value of a local operator after the quench 

(0(t, {r,})> = (y \e* m O({r t })e^ m \y ) (1) 

where |\l/o) is the initial state and H is the hamiltonian after the quench. The last 
relation can be written in path integral form and to assure its convergence one should 
first insert damping factors e~ eH with e — > as follows 

(0(t, {r t })) = Z- l ^ \e mt ^ H O{{n})e- lHt ^ H \m Q ) (2) 

where Z = (^ \e~ 2eH l^o) is a normalisation factor. If we analytically continue to 
imaginary time r then we obtain the same expression as that corresponding to a strip 
of width 2e where l^o) plays now the role of boundary conditions along both borders of 
the strip. Using the arguments of Renormalisation Group (RG) theory one can argue 
that as long as |\l/o) is translationally invariant, it can be safely replaced by another 
state that corresponds to an RG-invariant boundary condition, without changing the 
asymptotic behaviour of (2) in the limit e — > 0. In particular it turns out that for the 
bosonic field we discussed in the introduction and the deep quench limit, we have to 
impose Dirichlet boundary conditions (b.c.) forcing the field to vanish on the boundary. 
The same RG arguments give e a physical meaning as the typical time scale of the 
dynamics near the ground state of H , that is the inverse initial mass l/m . 
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This suggests that an inhomogeneous quench can still be formulated on a strip with 
Dirichlet b.c. but having variable width 2(e + h(x)) where h(x) expresses the variation 
of the initial mass 2(e + h(x)) ~ TOq (x). Now let us assume that the hamiltonian after 
the quench is massless, i.e. m = 0. In this case we can use CFT techniques and solve the 
problem directly for general h(x) by exploiting the conformal invariance of the theory. 
Indeed if we map the variable width strip (VWS) to the simpler geometry of a strip 
with constant width 2e (CWS) using a conformal mapping w — > z = g(w) then the 
transformation law of correlation functions of (primary) operators under such mappings 
will allow us to derive the corresponding expressions in the VWS from those in the CWS 
which are already known (Fig.l). The appropriate conformal mapping should satisfy 
the condition that 

lmg(x ± i(e + h(x))) = ±e (3) 




Figure 1. Schematic representation of the transformation from the VWS to the CWS. 



so that it deforms the boundaries as required. Notice that unlike the convention used 
in previous papers that the strip ranges from r = to r = 2e, we now center the 
strip in the middle so that r ranges from — e to +e. This is a well-posed mathematical 
problem whose solution can be found in general by solving Laplace's equation. In this 
way we obtain g{w) as an integral along the boundary of the VWS, a form that depends 
implicitly on h(x). This form is both cumbersome and actually inappropriate for our 
purposes since in order to extract useful information about physical quantities we only 
need to consider special asymptotic limits. To this end we will use two approaches: 

• the limit in which the transformation is infinitesimal, i.e. mo(i) varies only slightly 
in comparison with some average value. In this limit, if we define the inverse 
transformation as 1 + / = p -1 , the boundary conditions for / are to first order in h 

lmf(x ± ie) = ±h(x) (4) 

The solution can be found by using the general solution to the Laplace equation or 
by manipulating the analyticity of /. We then find that, up to an irrelevant real 
additive constant, / is given by 

+0O 

f{z) = J ds h{s)F{z - s) (5) 

— oo 
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with kernel 

1 1 / , TTZ \ 

F(z) = ; = — tanh — + 1 6 

v ; e(e-We + 1) 2e V It J K J 

where we have assumed without loss of generality that h(— oo) = 0. In Appendix 
A we explain all of the above in detail and present an elegant derivation of the 
solution. 

As we will see, to calculate the physical quantities we want it is sufficient to know 
the form of f(z) for real values of z. In addition if the distribution h(x) changes 
only slowly in comparison with e then, assuming that the following expression is 
convergent, the asymptotic form of f(x) in the limit e — > that we are interested 
in, is 

1 X 

f(x) = - J ds h(s) (7) 

— oo 

since F(x — s) can be written in this limit as Q(x — s)/e. 

the asymptotic behaviour of the transformation for large \z\, which turns out to 
determine the behaviour of physical quantities for large times and separations. Since 
the initial distribution h(x) must be bounded, we can distinguish two important 
cases: either it tends to the same value as x — » ±oo or to different values for each 
limit. In what follows we will consider one characteristic example for each case, 
focusing on the second which is more interesting: a bump distribution and a step 
distribution, both localised at the origin. Having assumed that h(—oo) = 0, for 
Rez — > — oo the transformation becomes asymptotically equal to the identity, that 
is f(z) — > 0. If we call h(+oo) = a then for Rez — > +oo the transformation should 
rescale the strip width by (1 + a/e) and since it is conformal must also rescale the 
x-direction by the same amount (Fig. 2). This means that the asymptotic form of 
f(z) as Rez — > +oo can only be 

f(z)~*z + (3 (8) 

where (3 is real so that, if a — 0, f(z) corresponds to a translation along the in- 
direction for Rez — > +oo. Overall the asymptotic form of the transformation from 
the CWS to the VWS is 

!z for Re(z) — > — oo 

(9) 
(l + a/e)z + P forRe(^)^+oo 

and that of the inverse transformation g 

| w for RefuO — > — oo 

g(w) ~ I K J 10 

\(w- P)/(l + a/e) for Re(iu) -> +oo 

The bump distribution corresponds to a = 0, while the step distribution to a ^ 
where a is the height of the step. It can be verified that the infinitesimal form of the 
transformation obtained from (5) and (6) has the expected asymptotic behaviour in 
these two cases. Notice that due to the fact that h(x) is supposed to be sufficiently 
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Figure 2. Conformal map from the CWS to a stepped width strip, as illustrated 
by the deformation of gridlines of the CWS. The transformation used is f(z) = 
aAlog (l + e 2z / A ) /2e where a controls the increase in the width in the right half and 
A the distance over which the change happens. Notice that due to the fact that the 
transformation is conformal, it scales the strip in the longitudinal direction by the 
same amount as in the transverse direction, in the right half. 



smooth in (5), we cannot use the Heaviside step function Q(x) to model the step 
distribution, since otherwise the result is incorrect. Also note that from (7) 

(3 = ~ Jim^ (J ds h(s) - xh(+oo) j (11) 

We can now calculate the evolution of physical quantities. We are mainly interested 
in the energy flow which is the off-diagonal component of the stress-energy tensor T 01 , 
the correlation function of (primary) field operators and the entanglement entropy. The 
results obtained for the strip should be analytically continued to real time r — > it and 
should give the correct asymptotic behaviour for times large in comparison with Mq^x) 
after the quench. 

In the next subsections we derive general formulae using the first approach. The 
derivation can be outlined as follows: we use the CFT transformation laws to express 
observables in the VWS in terms of those in the CWS in the form 



O vws ({ Wi }) = F [{g{wi),g(wi)}, O cws ({g( Wl )})\ (12) 
or for infinitesimal transformations, if we keep only the first order in / terms 



O vws ({wi}) = O cws ({wi}) + f(wi)}\VO cws ({wi}) (13) 

where the functionals T* and T involve derivatives of the functions g and / respectively, 
as well as their complex conjugates, with T being linear in /. Also note that w% = Xi+iTi 
and P is a suitable differential operator. 
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By the symmetry of the VWS under z ^ z we deduce that 

g(w) = g(w) and similarly f(w) = f(w) (14) 

which is satisfied by the kernel (6). This means that after the analytic continuation to 
real times r — > it, (13) gives for the observables in the inhomogeneous problem 

0({ Xi , U}) = O ({xi, U}) + F[{f{xi - U), f( Xi + t t )}]vO ({x t , U}) (15) 

where the subscript stands for the homogeneous expressions. This proves that it is 
only the restriction of / along the real axis that we need in our calculation. As a last 
step we take the limit e — > 0, also using (7) to express / in terms of h. 

2.1. Energy flow 

Let us first calculate the energy flow. We will do this calculation in detail in order 
to demonstrate the procedure outlined above. Assuming that f(z) is infinitesimal, the 
inverse map is to first order w — >■ z = w — f(w) and the well known CFT formula for 
the transformation of the stress-energy tensor in complex coordinates gives 

(T vws (w)) = (z'(w)) 2 {T cws ) - — {z,w} 

= (l-2f'(w))(T cws ) + y"'(w)) (16) 



where 



C / 7T x 2 



(T cws ) = (T cws ) = - {-) (17) 



The energy flow is given by the expectation value of T 01 (w,w) = i(T(w) — T(w))/2tt 
where we have taken into account the CFT normalisation factor — l/2n. Using the 
property (14) that is satisfied by / we find that 

= (1 " 2f(w))(T cws ) + y"\w)) (18) 
and combining the above results 

/ 7T \ 



£ (-/'H + /») + /'»- /'» 



(19) 



_V2e/ 

Now we substitute w = x + ir and analytically continue to real times r — > it to obtain 



(T°\x,t)) = 



2An 



{£) {-fix -t) + f'(x + t)) + f"'(x -t)- f'"(x + t) 



(20) 



Notice that, since T^ v is not a scalar but a tensor, T 01 (x,t) in real time equals 
—iT 01 (x,it) in imaginary time, that is apart from the substitution r = it we also have 
to multiply by dt/dr = —i. Finally taking the limit e — > and assuming that h(x) 
varies slowly in comparison with e so that we can use (7), we end up with 

f'TT 

( T01 M) = ^( h (x + t)-h(x-t)) (21) 

From the last equation, it is apparent that the energy flows from the initial time 
hypersurface to both directions in a wave-like fashion with speed equal to 1 (the 
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characteristic speed of the system, typically the speed of sound in condensed matter 
systems). The fact that it satisfies the wave equation is a general property of the stress- 
energy tensor in massless 2d theories, i.e. in CFT. What is non-trivial is the specific 
dependence on the initial distribution h(x). 



2.2. Correlation functions 

The next physically interesting quantity we can derive by the CFT transformation laws 
is the correlation function of a scalar primary field operator <3>. If we define 

C(z 1 ,z 2 ) = mz 1 )^(z 2 )) (22) 

then its transformation law is 

C vws (w 1 ,w 2 ) = \w'(z(w 1 ))w'(z(w 2 ))\~ x C cws (z(w 1 ),z(w 2 )) (23) 

where Wj = Xj + irj and x is the scaling dimension of $. Following the procedure 
described earlier we can write this transformation law in infinitesimal form. In particular 
to first order in / we have \w'(z(w))\ = 1 + Ref'(w). The CWS correlation function 
C CW s(zi, z 2 ) has been found in [8], but we can make a few more steps without using 
its explicit form, taking into account only the fact that it is invariant under space 
translations and interchange of the imaginary time variables. After some algebra we 
find that the equal time VWS correlation function to first order in / is given by 

d 



C vws (xi + ir,x 2 + ir) = ^1 - Re[/(xi + ir) - f(x 2 + ir)] 
-hm[f( Xl + ir) + f(x 2 + ir)}^- 



d(xi - x 2 ) 



-xRe[/'(xi + ir) + f'(x 2 + ir)]) C cws ( Xl + ir, x 2 + ir) (24) 
which after the analytic continuation to real times becomes 



C(xi, x 2 , t) = C (xi - x 2 , t) - 

_ l dC ( Xl -x 2 ,t) 
2 d(x 1 - x 2 ) 



(fix, -t) + f( Xl +t)- f(x 2 -t)- f(x 2 + t)) 



+ — (f(xi - t) - /(X! +t) + f(x 2 -t)~ f(x 2 + t)) - 

- ^C ( Xl - x 2 , t) (fix, -t) + f( Xl +t) + f'(x 2 -t) + f'(x 2 + t)) (25) 

where the subscript denotes the homogeneous case. Now we can use the explicit form 
of the homogeneous correlation function C (r,t) 

C(rt)= \(-Y cosh (W 2e ) + cosh(7rt/6) 1 x £ ^ f if t < r/2, 

olr ' ' \\2e) 8sinh 2 (7rr/4e)cosh 2 (7rt/2e)J ~ \ e ~x*r/2e tft>r/2. 

take the limit e — > and use (7) to finally find 

C(x 1 ,x 2 ,t) = C (x 1 -x 2 ,t) \l + ^ (Q(t- \x! -x 2 \/2) r * + r i+ * ds h(s)+ 

4e z \ Jx 2 -t Jx 2 +t 
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rxi+t rX2+t \ 

+e(|a;i -x 2 \/2-t) / +/ dsh(s))- 

J X\— t J X2—t / 

(h(x! -t) + h(x! +t) + h(x 2 -t) + h(x 2 + t)) 



(27) 



where we assumed that x\ > x 2 without loss of generality. Although the physical 
significance of this relation will become transparent later when we apply it to concrete 
examples, it is already clear that the effect of the inhomogeneity of the initial state 
propagates in waves with unit speed and that this gives rise to different behaviour 
inside and outside the horizon at t = \xi — x 2 \/2. 

2.3. Entanglement entropy 

We now turn our attention to the evolution of the entanglement entropy whose 
calculation turns out to be an application of the results of the previous subsection, 
since it can be expressed in terms of a correlation function of primary operators too. 
The entanglement entropy between a subsystem A defined by an interval [xi,x 2 ] of 
length / = \x\ — x 2 \ and the rest of the system, is equal to 

d_ 
dn 



S A = —TiaPa log Pa = - T^Trp^ 



(28) 



n=l 

where pa is the reduced density matrix of the subsystem A. From earlier work [25, 27] 
we know that Trp^ turns out to be proportional to the correlation function of primary 
field operators $±„ defined on the strip geometry and having complex scaling dimensions 
A n = A n 

(I _ - 1 1 _ _ I \ 2nA « 

1 1- 1 n 2 ~ n ^ ( 29 ) 
\ z l ~ z 2\\ z l — z 2\\ z l — z l\\ z 2 ~ z 2 1 / 

where z±, z 2 correspond to the edges x±, x 2 of the interval and imaginary time r and 

^iii'-h) (30) 

Notice that, since Ai = and Trp^ = 1, c\ must be equal to 1. It should be mentioned 
that although (29) is supposed to be valid for positive integer n the analyticity of the 
expression allows us to calculate its derivative with respect to n that appears in (28). 
From (28), (29) and the transformation law of correlation functions under conformal 
mappings z — > w(z) which is 

{QnM^-nM) = k , (^K))w , (^(w 2 ))|- 2nA "(<I»„(^K))<I»-„(^K))) (31) 

one can easily draw a general formula for the entanglement entropy in the inhomogeneous 
problem. Indeed, if we define the correlation function 

C(z u z 2 ) ee <* B (* 1 )<M*)> 1/(anA " ) = | u- 1 ~ "f 2 ~ =1 ( 32 ) 

\Z\ — Z 2 1 1 Z\ — Z 2 \\Z\ — Z\\\Z 2 — z 2 1 

then since C does not depend on n, the entanglement entropy is 

Sa = ~y (c n C 2nA ") = log C + const. (33) 

71=1 
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On the other hand, from (31) we see that C itself transforms as 

C vws (w 1 ,w 2 ) = \w\z(w 1 ))w\z(w 2 ))\~ 1 C cws (z(w 1 ),z(w 2 )) (34) 

which is a special case of (23) for x — 1- Hence in the infinitesimal inhomogeneous 
case C(x±,x 2 ,t) is given by (27) with x = 1. Substituting into (33) and using the 
homogeneous form of the entropy already known from [27] 

7r\ 2 cosh(7r//2e) + cosh(7rt/e) 
,2e) 8sinh 2 (7r//4e) cosh 2 (7rt/2e)_ ~ 

if t < 1/2, 



S A0 (l,t) = --]og 



~ 3 loge+ { 



cut 

Ik 
cnl 

12e 



(35) 



if t > 1/2. 

we obtain the entanglement entropy after an inhomogeneous quench 
S(x 1 ,x 2 ,t) = S (\x 1 - x 2 \,t) - 



C7T 



24e 2 



Xi+t X2+t 



Q(\xi-x 2 \/2-t)l J + J dsh(s)) + 



XI— t X2 — t 
X2—t X2+t 



+Q(t-\x 1 -x 2 \/2)l J + J dsh(s) 



\X\—t Xl+t 



+ 



+ 



12e 



(h(x 1 -t) + h(x 1 + t) + h(x 2 -t) + h(x 2 + t)) 



(36) 



Similar comments like those for the correlation function apply here. The wave- 
like propagation of the entanglement and the different behaviour inside and outside the 
horizon will be fully explained later using the quasiparticle picture. 



3. Application to the bump and step initial distributions 

Having found the general formulae we are now going to understand their characteristics. 
First of all, it is clear from (21), (27) and (36) that all quantities can be written as 
functions of the form f(x + t)±i(x — t). This means that the inhomogeneity of the initial 
distribution evolves in the form of waves to both directions, verifying the physical picture 
of quasiparticles emerging from the initial time hypersurface. Any initial inhomogeneity 
located for example at the origin, will only affect a distant point x after time t — \x\. 
This is another manifestation of the horizon effect that we already mentioned and a 
direct consequence of causality. 

To make our results more transparent we apply them to two main cases, the bump 
and step distributions, emphasizing the latter. As models of those cases we should use 
smooth distributions localised for example at the origin and preferably such that f(z) 
can be calculated in closed form at least from (7). The parameter (3 which corresponds 
to the size of the bump in the first case and a which is the height of the step in the 
second, must be both small if we wish f(z) to be infinitesimal so that we can use (5). 
On the other hand, results that do not rely on this restriction can be obtained from the 
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asymptotic approach using (9) and (10). These will be valid away from the horizon lines 
t — \x\. The last approach is more useful for the qualitative descriptions of this section. 

Let us start with the energy flow and assume for the moment that / is infinitesimal. 
The meaning of (21) is almost obvious: the initial energy density distribution propagates 
as in the classical wave equation. Notice however that this relation is only valid for slowly 
varying distributions relative to e and that in general the horizon will be smoothed over 
a distance of order e in the CFT approach. 

To visualise the evolution we will describe what happens in the above two cases. In 
the case of the bump, the latter splits into two equal parts each of which moves to the 
two different directions. In the case of the step, the energy flow is non-zero only inside 
the horizon, that is for t > \x\, where it takes the constant value C7ra/96e 3 . Notice that 
if the strip width is larger on the right than on the left, which means that the opposite 
is true for the initial energy density, then after the quench the energy flows to the right 
as it should. Of course as we can check using the asymptotic approach, this behaviour 
is correct away from the horizon even for finite transformations. 

Next comes the correlation function. We recall [7, 8] that in the homogeneous case 
(26) the correlation function Co decays exponentially in time until t — r/2 and then 
saturates to a value that depends exponentially on r and we will use the infinitesimal 
corrections to C (27) along with the asymptotic form (9). Let us concentrate on the step 
distribution and consider two points x\ and X2 separated by some distance r = \x\ — a^l, 
both lying on either the left or the right half of space, far away from the origin. If the 
pair of points is on the left then from (27) and (9) we see that right after the quench, C is 
equal to and evolves exactly like C since all the corrections vanish. This is because the 
two points have not yet been affected by quasiparticles from the right half. This starts 
happening at t — min{|xi|, ja^l} and C is changing until time equal to r has passed. 
If now the pair of points is on the right half then right after the quench the correction 
terms in (27) give exactly the first order corrections due to the substitution of e in Co 
by e + a. That is C has the homogeneous form Co but with the local strip width e + a 
as expected since the two points are not affected by quasiparticles from the left half. 
As before this happens at t — min{|xi|, \x2\} and for time equal to r. After this time 
C takes a value that is the same for the left and right half and equal to the average 
of the previous local values. To summarise, C evolves initially like the homogeneous 
correlation function Co corresponding to the local strip width until it enters the horizon, 
when the contribution of quasiparticles from both halves mixes to the average of the 
left and right value. 

A similar analysis holds for the entanglement entropy iSU(i). In the homogeneous 
case (35) this increases linearly with time until t = 1/2 when it saturates to a value 
proportional to the length I of A. Notice that there is also some constant amount of 
entropy that depends only on e in the form |cloge. This offset is the part that is already 
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there before the quench. Indeed since the theory is then massive the entropy according 
to [25] is |clog£ where £ is the correlation length, £ ~ Wq 1 ~ e. Moreover we mention 
that the entanglement in a massive theory is roughly speaking "located" close to the 
boundary points of A over a distance £. 

If we focus on the step distribution and consider again two cases for the position of 
the interval A, in the left or right half and away from the origin, then (36) tells us that 
Sa saturates first to the local homogeneous value at t = 1/2, then it starts changing at 
time t = min{|a;i|, \x 2 \} when the first quasiparticles from the opposite half enter A and 
finally after time equal to I it saturates again to its final value 

c , ncl ( a \ c a .„„. 

This value is the same in both cases and exactly equal to the average of the initial 
saturation values on the left and on the right. It is also the same even if A was in the 
middle containing the origin, as far as the length / is the same (Fig. 3). 

It is actually possible to go further and determine this asymptotic value of Sa for 
arbitrary finite h by using the general asymptotic form of the transformation g (10) 
together with (33) and (34). To this end we need however the expression for C cws 
corresponding to different imaginary times because the points x x + ir and x 2 + ir will 
be mapped by g to points with different in general imaginary parts T\, t 2 in the CWS. 
The required expression can be obtained from [8] 

C cws (r + iT!,iT 2 ) = 
it \ 2 cosh(7rr/2e) + cos(7r(r! + r 2 )/2e) 



2eJ 2(cos(7r(Ti - r 2 )/2e) + cos(tt(ti + r 2 )/2e))(cosh(7rr/2e) - cos(tt(ti - T 2 )/2t))^ 
When we analytically continue to real times and take into account (10) we find that in 



(39 a) 
(396) 
(39c) 




+ a/ej \ 1 + a/eJ 1 + a/e 
\xi — x 2 \- Finally substituting into (38), (33) and (34) we arrive at 

fe s ^'> = if \{\ + 7T^> + 1 (log£ + log(£ + a)) (40) 

This leads to the important conclusion that the entropy finally saturates to a value 
independent of the position. It only depends on the length of the interval /, no matter 
whether it is in the left or right part or somewhere in the middle. As far as the 
initial conditions are concerned, it is completely determined by the two limits of h(x) at 
±oo, being the average of the corresponding homogeneous values. Bumps or any other 
characteristics of h(x) do not affect the final saturation value. 
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Figure 3. Evolution of the entanglement entropy in the case of a step distribution, 
with the strip width smoothly varying over a distance of order e = 1, from 2 on the left 
of the origin to 4 on the right. The actual transformation used is f(z) = \ log(l + e 2z ) 
in which case the inverse transformation can be found analytically. The subsystem A 
under consideration is of length / = 10 and placed in 3 different positions with respect 
to the origin: with its middle at x m = —30 (i), (ii) and 30 (iii), all in units of 
e = l. We notice that in each case Sa first saturates to the homogeneous saturation 
value that corresponds to the local (i,iii) or average width (ii), at time t = 1/2 = 5. 
Note the different slopes and offset values before the saturation. In cases (i,m) and at 
time t = 25, i.e. the distance of the closest to the origin boundary of A, S4M starts 
changing again, since the first quasiparticles from the opposite half of space enter A. 
After time equal to / = 10, the entropy saturates again to its final value which is 
common for all positions and equal to the homogeneous value that corresponds to the 
average width. In the plot we used (33), (34) and (38) which are exact not only for 
infinitesimal transformations. 



Another interesting question is to find the entanglement entropy between the left 
and right halves in the case of a step initial distribution. Again we can do this for 
finite a. Now that the subsystem A is infinite, Trp^ transforms like ($ n (z)) where z 
corresponds to the edge of A (i.e. the origin) and time t. This is because the other edge 
is at infinity and does not contribute to the correlation function (<&„<!>_„). As before we 
can write Sa = — (c/6)log(7 where C(z) = ($ n (2)) 1 ^ 2nAn ^ transforms as 

CUH = \w'{z{w))\- l C cws (z{w)) (41) 
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and on the constant width strip takes the form found again in [8] 

6 - (x+iT) = h^k (42) 

where as always we have modified the result for a strip centered in the middle. Setting 
w = + ir in (41) and following the same procedure as before we finally find that for 
large times 

that is Sa increases linearly in time with a rate that is the average of the homogeneous 
rates corresponding to the two different strip widths. 



4. Free field theory approach 

Let us now study the inhomogeneous quench using a different approach: the real time 
evolution as prescribed by free field theory (FFT). This approach does not have the 
restriction that the hamiltonian after the quench be massless and it offers both a validity 
check of the CFT results in the massless case and an extension of the results to the 
massive case. As we will soon see, the time evolution of the two point correlation 
function can be easily found by solving the Heisenberg equations of motion. Then 
the problem reduces to finding the initial correlation function, that is the correlation 
function in a theory with spatially inhomogeneous mass. This is formally equivalent 
to a scattering problem where the mass plays the role of the potential and it can be 
generally solved by means of perturbation theory, assuming that the spatial variation 
rj(x) of the mass is small relative to some characteristic value m . We are interested in 
the first order correction in rj, which is sufficient in order to make a comparison with 
the conformal result. Unlike the CFT method where the correlation functions and the 
energy momentum tensor were derived independently, we will first calculate a general 
expression for the correlation function of the field and from this the energy flow. Then we 
explore the asymptotic behaviour of the latter in the special case of a step distribution 
in the deep quench limit. 

4-1. Perturbative solution and calculation of the propagator 

We start with the calculation of the two-point correlation function of the field operator 
<f). The Heisenberg equation of motion for <fi after the quench at t — is 

4>(x,t) = (%-m 2 )<f>(x,t) (44) 

which can be readily solved in Fourier space where it reads 

4>{k-t) = -u 2 k 4>{k-t) (45) 

with uj 2 = k 2 + m 2 and (f)(k;t) = J dx e~ lkx (f)(x , t) the Fourier transform of (f)(x,t). The 
solution is 

(j)(k; t) = (f)(k; 0) cos u k t + <j)(k; 0) Sm — (46) 



Inhomogeneous Quantum Quenches 



16 



so that the two-point correlation function is 

(0(xi,tiMa; 2 ,f 2 )) = J ^^e^ 1+ife2a;2 ((0(A; 1 ;O)0(A; 2 ;O))coscc; 1 t 1 cos^ 2 + 
+{<t>{k 1 - O)0(A; 2 ; 0)) cosc^ 811 ^ 2 + O)0(A; 2 ; O)) 8 "^ cos^ 2 t 2 + 

w 2 CJi 

/i/i nw/i rvw sinc^iti sinu; 2 t 2 \ 
+ (0(fci;O)0(fc 2 ;O)) (47) 



where Ui = u^- 

We have isolated the time evolution and we now need to calculate the correlation 
functions just before the quench, when the system lies on the initial ground state. The 
expectation values (...) are therefore meant to be evaluated on this state so that we can 
forget the quench for the moment and focus on the field theory before that. If we write 
the initial mass distribution as mo(x) = m + rj(x) and assume that \rj(x)\ <C m then 
the field equation before the quench is 

(02 _ d 2 + (mQ + rj^f)^ t ) = (48) 

or to first order in i](x) 

(d 2 -dl + m 2 + 2m r](x))(j)(x, t) = (49) 

For later comparison with the CFT results, note that the correspondence between the 
strip width and the initial mass leads to e = l/m and so h(x) = —r}(x)/ml where the 
variation is always assumed relatively small. 

The field equation above can be solved perturbatively in r](x). The method is same 
as that used in the derivation of the Lippmann-Schwinger equation in scattering theory, 
with the inhomogeneous terms corresponding to the potential. The homogeneous field 
equation is the well-known Klein-Gordon equation 

(d*-d 2 x + m 2 )M^t)=0 (50) 
whose retarded Green's function, satisfying 

(d 2 -d 2 x + m 2 )G 0R (x, x', t, t') = -i5(x - x')5(t - t') (51) 

is 

G 0R (x, x', t, t') = G(t - t')(G (x -x',t- t') - G (x' - x, t' - t)) (52) 

with 

Go{x _ x ' jt _ t >) = ^t)4> G {x' ,t')) = [ % ik ^ x ')"*-*(*-*') J- (53) 

J Ztt z&ok 

where cook = \Jk 2 + m\. From these last relations, (49) can be written in the integral 
form 

<f>{x, t) = <f) (x, t) - 2im J dt' J dx'G 0R (x, x', t, t')r](x')(j)(x' , t') (54) 
and since r\ (x) is small this leads to first order to the solution 

<j>(x,t) = (f) (x,t)-2im I dt' [ dx'Go R (x,x',t,t')r](x')(t) ( x '^')+ (v 2 )(^) 



Inhomogeneous Quantum Quenches 



17 



<£(* i; O)0(A; 2 ; 0)) = 2^^ + k 2 )^- + 2 ™«f + **) ( J_ _ M (58a) 



from which follows that before the quench 

(<l)(x 1 ,t 1 )(f)(x 2 ,t 2 )) = (4>o{xiM)<t>o{x2,h)) ~ 

-2im J dt' J dx'G OR (x 1 ,x',ti,t')r](x')((f) (x',t')(t)o(x2,t2)) - 

-2im J dt' J dx'G 0R {x2,x , ,t 2 ,t')rj(x'){(l) (xi,t 1 )(l) (x , ,t')) + O(ri 2 ) (56) 

From (52) and (53) we can work out the time integral (adding appropriate small 
imaginary numbers to the lower integration limits to ensure convergence). Passing 
to momentum space, the final result for the correlation function to first order in r\ is 

$(ki;tS(k 2 ;t 2 )) = 

e -^oi(ti-t 2 ) 2m fj(k 1 + k 2 ) / e -*"oi(*i-fc) e -^o2(ti-t 2 )\ 

= 2n6{k 1 + k 2 ) - + —r^ - - 57 

2u 01 (woi - ^02) V 2 ^oi 2^02 / 

where u 0i = cj fc; and ^(fc) is the Fourier transform of 77(2;). 

Going back to (47), we see that we need the following equal time correlation 
functions between field and conjugate momentum operators n = <p, which can all be 
found from (57) by applying appropriate time derivatives 

1 2m fj(k 1 + k 2 ) / _Jl 1 

jqi uj^ — luq 2 \2o;oi 2cjo2^ 

(fe; O)0(A; 2 ; 0)) = nid^h + k 2 ) (586) 

(^(fci; O)0(A; 2 ; 0)) = — 7ri<5(fei + k 2 ) (58c) 

(0(fci; O)0(/c 2 ; 0)) = 2vr5(A; 1 + A: 2 )— + 2 2 o~ ( 58rf ) 

Z C^Qi — Cc?Q2 \ Z Z / 

As a check, we can verify the consistency with the canonical commutation relations 
[4>{ki] 0), 4>{k 2 ] 0)] = 2m5(k 1 + k 2 ). The next step is to substitute into (47) to find that 
the propagator after the quench is 

(<l)(x 1 ,t 1 )(j)(x 2 ,t 2 )) = 

= ( — e lk ( x i~ x i) sincjti sincjt 2 + cos cut 1 coscut 2 — sincj(ti - t 2 )\ + 

J 2n \2uj z 2ujq 2uj J 

+mQ f dkl ^ e ifci3!i+ifc3x 3 g(^i + fc 2) / sin u;itisincj 2 t 2 _ cos^iticos u 2 t 2 \ 

i 27T 27T C^oi + ^02 V ^1^2 ^01<^02 / 

The first line is just the propagator in the homogeneous case found in [8] and the second 
line is the first order correction due to the inhomogeneity. In the deep quench limit 
m — > oo, (59) can be written as 

<«*..*.)#(*.*>» =J~ /^■■-'-'^-%(») Sb ^^ 1 (60) 

This expression can be found directly from (47) if we notice that the highest order in 
mo contribution comes from the initial momentum- momentum correlation function (58 d) 
which can now be written in real space as (<fi(xi, O)0(x 2 , 0)) ~ m (xi)5(xi — x 2 )/2. This 
is because the initial correlation functions fall exponentially fast over a distance of order 
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I/tuq — > and so they can be approximated by 5-functions with suitable coefficients 
which, in the current inhomogeneous problem, depend on the "local mass" mo(x). 

In the language of Feynman diagrams (60) says that the dominant part of the 
correlation function between two points (xi,ti) and (x 2 ,t 2 ) comes from a diagram with 
two lines connecting the points with another one (s, 0) at the initial time hypersurface. 
What such a diagram represents physically is a pair of quasiparticles emerging from 
this hypersurface and reaching the two points. In general in order to cause correlations, 
two quasiparticles should emerge not necessarily from the same point, but from points 
separated by a distance of order l/m (s), that is the "local correlation length", and the 
deep quench approximation consists in saying that these two points are so close that 
can be effectively identified. 

In the massless case i.e. when m = 0, Uk = \k\ and the Fourier transforms in (60) 
can be calculated explicitly giving the correlation function as a convolution of the initial 
mass distribution 

(<f){x 1 ,t 1 )<j)(x2,t2)) = ^Jds m (s)©(ti - |xi - s|)0(t a - \x 2 - s\) (61) 

Notice that in the homogeneous case m = const, and for equal times, we recover the 
deep quench massless propagator discussed in [8] 

fo if * < \x 2 -x 1 \/2, 

(0(x 1 ,t)0(x 2 ,t))= 17 ' (62) 

[ (2t - \x 2 - xi\)m /8 if t > \x 2 - Xi\/2. 

The meaning of (61) is quite transparent: the correlation function between two 
spacetime points depends on the values of m (s) only at those points of the initial 
time hypersurface that lie inside the horizons of both points (Fig. 4). The physical 




Figure 4. Physical explanation of (61). The connected two-point correlation function 
is determined by the values of the initial distribution at the overlap of their horizons 
(red thick line). 



interpretation of correlations will be fully explained later using the quasiparticle picture. 
Also it will be shown that (61) is in agreement with the CFT formula (27). 
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4-2. Energy flow 

The energy flow after the quench, i.e. the expectation value of the T 01 component of the 
stress-energy tensor, can be obtained by acting on the correlation function with —d tl d X2 
and then setting x 1 = x 2 and t\ = t 2 

(T 01 (x,t)) = -{n(x,t)d x <i>(x,t)) = - lim ]imJ--?-{cf>(x' ,t')<l>(x,t)) (63) 

x 1 — *x t'—>t at ox 

From (59) we then find 

T°\x,t) = -m i /^^(W*±M ta (- cosc^shu^ + sin cos u; 2 A 

J 2n 2n UJ 01 + UJ 02 \ U0\002 ^01^02 / 

The homogeneous term, as expected, does not contribute to the energy flow. We are 
mainly interested in the value of the energy flow at the origin x = and for t — > +oo. 
As before, in the massless case the expression above simplifies significantly 

mm , . r dki dko fiik] + ko) ( k\k 2 \ , s . . 

T 01 0,0 = ~m i / — \ 1 + sin (h + k 2 )t 65 

J 2tt 2vr 2(cjoi + cj 02 ) \ ^01^02/ 

but unlike before, we cannot use the deep quench limit in the first place as in (60) 
as this would lead to ultraviolet divergences. In fact this should be expected since the 
double differentiation of the propagator dropped down a factor of k 2 uii making T 01 more 
sensitive to ultraviolet divergences. Therefore the approximation that m k, q is not 
valid anymore. 

To extract the large t behaviour of T 01 (0,t) we can probe the small frequency 
behaviour of its Fourier transform T 01 (0;uj) = J dt e~ luJt T m (0, t) instead. Notice that 
the time integration runs from —00 to +00 although (65) is physically valid only for 
t > 0. Since (65) is an odd function of t, T 01 (0; uj) must be an odd function of uj. This 
means that if T 01 (0, t) tends to a non zero value as t — > +00 then T 01 (0; a;) must behave 
like i/uj for uj — > 0, while if it tends to zero then T 01 (0;cj) will have no singularity at 
uj = 0. Setting k\ = k + q and k 2 = k — q 

T*(0;u>) = -moi [ ^ ^ Jl+ *>* " Q - ) [ dte~^ sin 2kt 

J 2lX 2lX [uj^ k+q ) + ^O(k-q)) \ ^0(k+q)^0(k-q) J J 

f dk dq f)(2k) / k 2 - q 2 \ n 
= -mo / t^t^tt, r 1 + 2n(5{2k 

J 2-K 2-K 2{UjQ( k+q - ) +0J (k-q)) V <^0(k+q)<^0{k-q) J 

= -\rno(ij(u) - fj(-oj)) f ^-j- ^— r 

4 J Zn (UJ ^/2+q) + Wo( w /2-q) ) 

and for uj — > 




T o1 (0;uj) ~ --m (fj(u) - fj(-u)) / £^ = -— m (fj(uj) - fj(-uj)) (67) 

J 0q 

Applying our result to the case when the initial mass distribution is a step function of 
step a, where fj(k) = ai/k we find 

T°i( ;^0)~-^ (68) 
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which, according to our discussion above, means that the energy flow at the origin tends 
for large times to a non-zero value, more specifically to — am Q /(8n). This result proves 
once again the non-diffusive behaviour of the energy flow. Notice that if a > 0, i.e. 
the initial mass and consequently the initial energy density is higher in the right half 
of space than it is in the left, then the minus sign indicates that the energy flows from 
right to left, as it is supposed to do. 

The massive case is also interesting independently from the massless one since 
thermalisation occurs for a different reason in the two cases. Following a similar 
calculation presented in Appendix B we conclude again that the energy flow at t — > +00 
tends to a non-zero value. 



5. Comparison between the conformal and free field theory results 



Since we were able to calculate the energy flow and correlation function both in the CFT 
and FFT approach, we can now compare our results. We start with the energy flow. In 
both cases we have shown that it does not decrease with time inside the horizon. If we 
compare the asymptotic values at the origin for large times using the correspondence 
relations e = l/m and h(x) = —7]{x)/rn^ we notice that the two results agree up to a 
numerical factor of 7r 2 /12. Consequently the CFT method is consistent with FFT as far 
as the qualitative behaviour of the energy flow is concerned, the only difference being in 
numerical factors. This remark holds in the homogeneous case as well if one calculates 
the energy density instead, as explained in Appendix C. 

Now for the correlation function we have to take into account that in the gaussian 
model (free boson) the primary field $, whose correlation function C(xi,x 2 ,t) is given 
by the CFT expression (27), is not the field of FFT but its imaginary exponential. 
More specifically C(x\,x 2 ,t) should be compared to ^e ig ^ xl ' t ^e~ iq ^ X2 ' t ^) where q is an 
arbitrary constant. Using properties of gaussian integrals we have 



(e 



iq<p(xi,t) e -iqcf>(x2,t)\ 



= e 



■%-<(*(*i.*)-*0w,t)) a > 



= exp 



-(G(x!,xi, t) + G(x 2 , x 2 , t) - 2G(xi, x 2 , t)) 



(69) 



where G(xi,x 2 ,t) is the correlation function of in the massless case, as given by (61). 
After some algebra we find that to first order in rj(x) 



(e 



iq<j>(x!,t) -iq<t>(x 2 ,t)\ 



o q 2 {G (x 1 -x 2 ,t)-G (S),t)) 



2 ~" 

" 16 / dS V{S) 



x 



X 



(0(t - \xx - s\) + G(t - |a:2 - s\) - 29(t - \xi - s\)0(t - \x 2 - s|)) 



(70) 



where G (xi — x 2 ,t) is the homogeneous correlation function (62). 

On the other hand, the CFT result (27) can also be written in a comparable form 



C(x u x 2 ,t) = C (x 1 - x 2 , t) 



1 + 



X7T 

4e 2 



J ds h(s)K(xi,X2,t; s) 



(71) 
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where the kernel K{x\,x 2 ,t\ s) is 

K(x u x 2 , t; s) = e(|xi - x 2 \/2 - t) (Q(t - \ Xl - s\) + S(t - \x 2 - s\)) + 

+Q(t - \x! - x 2 \/2)(e(s -x 2 + t)e(-s + x!-t) + e(s -x 2 - t)e(-s + Xl + 1)) (72) 

and we kept only first order in e terms. The homogeneous correlation functions 
exp [q 2 (Go(xi — x 2 , i) — Go(0, t))] and Cq{x\ — x 2 , t) have been already shown to be equal 
[8] assuming that q 2 = 4xvr. What we wish to verify now is that the first order corrections 
agree too. This is true since the integration kernels in (70) and (71) are in fact identical 
and even the numerical coefficients are equal if we use the previous substitution for q. 

6. Discussion and conclusions. 

6.1. The quasiparticle interpretation 

Our results can be easily interpreted using the physical picture of quasiparticles 
developed in earlier work [8] and verified for the special case of a domain wall initial 
distribution in [22]. This can be seen in the formula (36) that gives the evolution of the 
entanglement entropy for example. 

First we repeat how this physical picture applies to the homogeneous case. As 
always we consider a subsystem A which is an interval of length / while the complement 
is the subsystem B. Then the entanglement entropy between the two subsystems after 
the quench is given by (35) up to an additive constant independent of e. The constant 
value |cloge is the part that corresponds to the massive field theory before the quench. 
In that case correlations are restricted to between points separated by distances of order 
I/tuq ~ e which is the correlation length. Hence the entanglement between A and B 
is only due to correlations between points close to the two boundaries of A. Since the 
initial state has energy density much higher than the ground state of the Hamiltonian 
after the quench, it acts as a source of quasiparticles moving with unit speed to both 
directions. Two quasiparticles emitted from the same initial point are entangled and if 
they move to different directions they can reach the two different subsystems increasing 
the entanglement between them. The latter must be at any time proportional to the 
number of such entangled pairs which is simply proportional to the length of the region 
€(A) where these pairs are emitted from (Fig. 5). This is the total length of the left and 
right projections of A on the initial line apart from their intersection, since quasiparticles 
coming from there both end up inside A. Also the number of quasiparticle pairs emitted 
from some point must be proportional to the initial energy density at this point, which 
in the homogeneous case is constant and proportional to 1/e. 

Keeping all these in mind, we can easily see that shortly after the quench, the pairs 
that contribute to the entanglement are emitted from within distance equal to t from 
the boundaries of A. This causes a linear increase ~ 4t/e. However at time t = 1/2, the 
region S(A) reaches its maximum length 21, since the left and right projections of A on 
the initial line no longer overlap. This explains the saturation to a value ~ 21 /e. By 
comparison with (35) the proportionality factor turns out to be equal to c7r/24. 
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Figure 5. Illustration of the physical interpretation of the entanglement entropy 
evolution using the concept of entangled quasiparticles. 

(i) Pairs of entangled quasiparticles emitted from the same point on the t = 
hypersurface. One of the quasiparticles of the pair denoted by dark red colour 
is inside subsystem A at some time t, while the other is in its complement B. 
Therefore this pair contributes to the entanglement between A and B 
at that time. In contrast, pairs like those denoted by light grey colour whose 
quasiparticles arc both in the same subsystem, cither A or B, do not contribute 
to S A (t). 

(ii) The red thick lines denote the regions on the t = hypersurface where entangled 
quasiparticles that contribute to Sa (i) come from. For t < 1/2 these regions are 
[x\ —t,x\ + t] and [x2 — t,X2 + t], while for t > 1/2 they are [ari — t, X2 — t] and 
[xi +t,x 2 + i\. 
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Following the same arguments it should be straightforward to generalise to the 
inhomogeneous problem. From (36) the part of the entanglement entropy that 
corresponds to the massive theory before the quench is 

S( Xl ,x 2 ,0) = Uoge + ^-{h{ Xl ) + h{x 2 )) (73) 

This is simply the sum of the contributions of the two boundaries, if we take into account 
that now it is the local correlation length at each point that should be used 

°- [log(e + h(x 1 )) + log(e + h{x 2 ))] (74) 
6 

and that h(x) <C e. To calculate the contribution of this initial entanglement after the 
quench, we must take into account that due to the wave propagation, A is affected not 
by the initial value at xi,x 2 but at their projections x\ ±t,x 2 ± t. This justifies the 
term 

-^-(h(x! -t) + h(x! + t) + h(x 2 -t) + h(x 2 + *)) (75) 

Now we calculate the contribution of the quasiparticles. As before, until t = 1/2, this 
comes from the intervals [x\ — t,x\ + 1] and [x 2 — t,x 2 + 1], but now the emission rate 
at some point s is proportional to the local initial energy density ~ l/( e + M s ))- This 
leads to 

Xl+t _^ x 2 +t \ 

I ds 7TW) + I ds ^W)) (76) 

\Xl-t X 2 -t K 'I 

and since h(x) <^ e 

Xl+t X2+t \ 

dsh(s)+ f dsh(s)\ (77) 



cut or 



CTll C7T 

12e ~ 24? 



6e 24e 2 

\xi— t X2~t 

The first term can be recognised as the homogeneous part while the rest is due to the 
inhomogeneity. Similarly after t = 1/2 we find 

X2~t X2+t \ 

J dsh(s)+ J dsh(s)) (78) 

KXl—t Xl+t / 

Adding altogether we obtain (36). 

A similar analysis applies to the expression for the correlation function C(xi, x 2 ,t). 
This shows that the quasiparticle interpretation successfully explains all the details of 
the evolution of correlations following an inhomogeneous quantum quench. This picture 
should also be valid in the massive case, the difference being that the quasiparticles now 
have finite lifetimes and propagate with various velocities up to the maximum one [8]. 
This smoothes out the horizon and leads to spatial oscillations of correlations inside and 
exponential decay outside it, as can be shown by stationary phase arguments. 
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6.2. Comparison to the thermal analogue 

Our findings enable us to answer the main question that motivated the study of this 
problem, that is if there are any similarities between an inhomogeneous quantum quench 
and its thermal analogue. This is the evolution of an inhomogeneous initial temperature 
distribution uq{x) which is given by the heat equation d t u(x,t) = Kd^u(x,t) with initial 
condition u(x, 0) = uo(x). The solution is 

+oo 

u(x, t) = . / dsuo(s) exp [— (x — s) 2 /4/tt] (79) 

— oo 

from which we can find the heat flow by Fourier's law 

dO A + °° 

j q (x,t) = ^r = -Xd x u(x,t) = - -j= J dsu' (s)exp[-(x- sf/Ant] (80) 

— oo 

As can be seen from the last equation and mentioned in the introduction, for any step- 
like initial distribution the heat flow at the origin decreases as 1/ \fi for large times 

Aa 

^"-va (81) 

where a is the size of the step. This can be considered as a characteristic of diffusive 
behaviour thus providing a simple test for our analysis. 

In a quantum quench on the other hand, the energy density which is supposed 
to be proportional to the effective temperature, if we could assign a local meaning to 
it, exhibits wave-like non-diffusive behaviour instead as we can see from (21). The 
same message comes from (68) and (B.6) when we see them through the lens of the 
aforementioned test. We thus conclude that the effective temperature is not meaningful 
as a local quantity. Notice that the non-decreasing of the energy flow is also true in 
the massive case m ^ 0, even though the quasiparticles have finite lifetimes. Recall 
that for m ^ it is the bosonic propagator itself that thermalises, unlike the conformal 
case where thermalisation occurs on the level of correlation functions of primary field 
operators. Thus the nature of thermalisation is qualitatively different in each of these 
two cases and one should consider it as an independent effect. 

Another aspect of the comparison is whether the entanglement entropy S en t 
resembles the thermodynamic one S t h- We already know from the homogeneous case 
that S en t becomes extensive when it saturates and of course the same holds in the present 
inhomogeneous case. One might be tempted then to define an entanglement entropy 
current and ask if this plays a role similar to the heat current j q = T{dS t h/dt) where 
T is the temperature. Once again this question is simplified for a step distribution and 
with the complementary subsystems A and B between which S ent is measured, being 
the two halves of space, on the left and right of the origin. From (43) we see that in CFT 
the entanglement entropy rate for large times is constant, not decreasing as it should 
happen if it exhibited diffusive behaviour like the heat current. 

As a final remark we will discuss the second interesting question arisen in the 
introduction which was what information about the initial state survives in the 
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stationary values of local observables for large times. As must be clear from the CFT 
expressions (21), (27) and (36) in the limit t — > +oo the only relevant parameters of 
the initial distribution h(x) are its asymptotic values in the limits x — * ±oo. This must 
be obvious considering the wave-like nature of the evolution. We should mention by 
the way that, even if the evolution was of diffusive nature, those two limits would still 
be the only relevant parameters determining the large time behaviour. For example 
the uniform relaxation temperature is simply the average of the initial temperature 
distribution / ds u (s) which equals the average of the above two limits «o(±oo). 
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Appendix A. Derivation of the conformal map 

We seek a conformal transformation g(w) from the strip with variable width 2(e + h(x)) 
to the strip with constant width 2e, which means that g must satisfy the boundary 
conditions (3). As mentioned in the main part of the article, these conditions have the 
disadvantage of being defined along the boundary of the VWS whose shape is nontrivial. 
Let us define the inverse transformation 1 + / = g^ 1 from the CWS to the VWS. If 
we expressed the problem in terms of / then the boundary conditions would be defined 
along straight lines, but (3) cannot be translated into sufficient boundary conditions for 
/ unless we assume that the transformation is infinitesimal. More specifically, if h(x) 
is small compared to e and sufficiently smooth, then g is close to the identity and / is 
of the order of h so that the boundary condition implied for / to first order in h is (4) 
lmf(x ± ie) = ±h(x). 

Now the problem reduces to finding the analytic function f(z) defined on the CWS 
and obeying (4). Imf is a solution of the Laplace equation which ensures that the 
boundary conditions above along with conditions at ±oo are sufficient for determining 
Im/, while Re/ can be determined using the Cauchy-Riemann equations up to a real 
additive constant. This is irrelevant for our purposes since we can check that it would 
vanish from all of the expressions for the physical observables that we are interested in. 
Hence we can fix its value by requiring that / is equal to zero for vanishing h. Without 
loss of generality we suppose that h(—oo) = and then 

lim f(z)=0 (A.l) 

Hez^—oo 

One way to find f(z) is using a generalised form of Cauchy's integral formula with 
an appropriate kernel. In particular we have 




(A.2) 
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where C is a closed contour around z and D(z' — z) is any function with a simple pole 
with unit residue at z' — z as the only singularity within C . If we choose C to go along 
the boundaries of the strip (Fig.Al) then 



+ie +h(x) 



00 


z C , 

• 


+ 00 


1 







-ie -h(x) 

Figure Al. Contour integration for the derivation of f(z). 



f(z) = - — : [ (D(s — ie — irj — z)f(s — ie) — D(s + ie + irj — z)f(s + ie)) ds + 

Z711 J 

— oo 

i 7 

+ — / (D(+oo + is- z)f (+oo + is) - D(-oo + is - z)f (-00 + is)) ds (A. 3) 

— e 

where rj — > 0+ simply serves as a reminder of the prescription that the pole is always 
enclosed inside the contour, even for z on the boundary. Using (A.l) and the symmetry 
of the problem under reflections z — » z which implies that f(z) = f(z), we find 

+oo 

f(z) = - — ; / (D(s — ie — irj — z) — D(s + ie + irj — z)) Re/(s + ie)ds — 
2m J 



— oo 

+oo 



— — f (D(s — ie — irj — z) + D(s + ie + irj — z)) lmf(s + ie)ds + 

Z7T J 

— oo 

i 7 

+ — / D(+oo +is- z)f(+oo + is)ds (A.4) 

— e 

The last relation suggests that we choose D(z) such that D(z + 2ie) = D(z) and 
\im Rez ^ +OQ D(z) = 0, because then f(z) is completely determined by the boundary 
condition (4) 

+oo 

f(z) = / D(s + ie - z)h(s)ds (A.5) 

71 J 



— oo 



Now it is relatively easy to see that one function that satisfies all the above conditions 
is 

D(z) = -7 j- — —r (A.6) 
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+00 

/(*) = " / e -« z -, )/e + 1 K°)d* (A.7) 

—00 

i.e. (5) and (6). 

Appendix B. Free field theory calculation of the energy flow in the massive 
case 



We will calculate the large time asymptotic behaviour of T 01 (x,t) at x = and 
with m^O following the same method as for the massless case, that is using its Fourier 
transform with respect to time T 01 (0;u). As before, if this is proportional to i/u for 
uj — > then T 01 (0, t) tends to a non zero value as t — > +00. Starting from (64) we have 



T 01 0,0 = -m i / tt^-^ o, , \ fc 2 o;i + ) sm^ + u 2 )t+ 

J 2n 2-k 2{uj 01 + uj 02 ) LW iu;o2 cjicj 2 / 

+ ( — — ) sin(cji - cu 2 )t 



(B.l) 



x 



Its Fourier transform for small frequencies uj is 
T (0;w) — / dk x dk 2 - ■ rk 2 u>i 

o7T J (^01+^02) W iCJ 2 0J X 0J 2 , 

X ^5(c<Ji — UJ 2 — Uj) — 5(oj\ — UJ 2 + Uj)^j (B.2) 

Notice that we have skipped the two ^-functions of the form 8{uj\ + uj 2 ± uj) that come 
from sin(c<Ji + uj 2 )t since they vanish for uj — > as uj\ + uj 2 > 2m > 0. Let us specialise 
in the step distribution 

, x f-a/ 2 if a; < 0, . . 

" W i + a/2 if*>0. (B3) 

whose Fourier transform is ^(/c) = ai/k. If we symmetrise the integrand under k\ — > — fci, 
k 2 —*—k 2 and &i <-> fc 2 we finally obtain 

T ul (0;w) ~ — — / dk x dk 2 - ■ -— 2 )b{uJi- uj 2 - uj) (B.4) 



57T J (CJ 1 + CJ 2)(^l - kf) \UJqiUJ 02 UJ\UJ 2 , 

Now the 5-function prescribes that uj\ = uj 2 + uj and substituting to the integrand we 
find to first order in uj 

f>;^0)~- m ° a * (m °- m2) fdk 1 dk 2 ^;f ) 5(uj 1 -uj 2 -uj) (B.5) 

16tvuj J ^01^1 

This shows that T 01 (0; uj) is of order 1/uj for close to zero, which means that as before 

T 01 (0,t) does not tend to zero for large times. Evaluating the integrals carefully we end 

up with the simple relation 

f>;c^0)~- (m °- m)a * (B.6) 

Attuj 

For m = we recover the already found result (68). Notice that the expression above 
gives zero if m — m . 
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Appendix C. Comparison of the CFT and FFT results for the energy 
density after a quantum quench 

We will calculate the energy density after a homogeneous quench of the mass from m to 
zero, using the CFT and FFT methods. The CFT result is obtained by the expectation 
value of the T 00 component of the stress-energy tensor on a strip of width 2e 

(TZ P ) = -^^2 ( ai ) 

The analytic continuation from imaginary to real time reverses the sign and when we 
also set c = 1 and take into account the correspondence relation e = l/m , we obtain 

(T°°) = ^ (C2) 

In free field theory the energy density (T 00 ) can be obtained from the propagator 
in the same way as the energy flow (T 01 ), but here we will show a slightly different 
derivation. The energy of the system which is conserved, is obviously equal to the 
expectation value of the hamiltonian after the quench H = | / (d<p) 2 dx with respect to 
the ground state \^f ) of the hamiltonian before the quench H = \ J ((<90) 2 + ml<f) 2 )dx. 
If we decompose the latter as follows 

(*o|#|*o> = (*o|(# " #o)|*o) + (^o|#o|^o> (C3) 
then the first part is 

- l 2 ml J dx (^ o |0 2 |^o) (C.4) 
which is easy to calculate since 

<w»> = /Sft^ (c ' 5) 

is a loop of the Feynman propagator with mass m in (l+l)-d, while the second part 
is the reduced free energy of a system with hamiltonian Hq = H + \m\ J <p 2 dx where 
the last term is considered as a perturbation over H . The reduced free energy per unit 
length in this case is 

f{ml) = llwf I log {e + ml) - log fc2 l (c - 6) 

Putting all these together we find that the energy density is equal to 



d 2 k 



log (1 + ml/k 2 ) - 



rrir 



(C.7) 



k 2 + ml 

The last integral is both ultraviolet and infrared convergent and if we set k 2 = um\ and 
integrate by parts, it gives 



m\ r du ml 



8tt J (1 + u) 2 8tt 

o v ' 

r 2 



(C.8) 



Comparing the two results we notice that they differ by a factor n /12, exactly as 
for the energy flow. 
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